Velocity of vortices in inhomogeneous Bose-Einstein condensates 
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We derive, from the Gross-Pitaevskii equation, an exact expression for the velocity of any vortex 
in a Bose-Einstein condensate, in equilibrium or not, in terms of the condensate wave function at 
the center of the vortex. In general, the vortex velocity is a sum of the local superfluid velocity, 
plus a correction related to the density gradient near the vortex. A consequence is that in rapidly 
rotating harmonically trapped Bose-Einstein condensates, unlike in the usual situation in slowly 
rotating condensates and in hydrodynamics, vortices do not move with the local fluid velocity. We 
indicate how Kelvin's conservation of circulation theorem is compatible with the velocity of the 
vortex center being different from the local fluid velocity. Finally we derive an exact wave function 
for a single vortex near the rotation axis in a weakly interacting system, from which we derive the 
vortex precession rate. 

PACS numbers: 03.75.Hh, 05.30.Jp, 67.40.Vs, 67.40.Db 
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INTRODUCTION 



Vortices in rapidly rotating trapped Bose-Einstein con- 
densates [H 13- SOI do not move with the local fluid flow, 
despite Lord Kelvin's circulation theorem which implies 
that a vortex in a perfect fluid generally cannot escape 
from a contour co- moving with the fluid 0, . To un- 
derstand how vortices do in fact move in rapidly rotating 
condensates, we derive here an exact expression for the 
velocity of a vortex, applicable both to equilibrium and 
non-equilibrium situations involving arbitrary numbers 
of vortices, through analyzing the motion of the vor- 
tex singularity directly from the Gross-Pitaevskii (GP) 
equation (Eq. J2J) below). This differential equation, a 
generalization of the Schrodingcr equation, determines 
the time dependence of the order parameter 'J (the con- 
densate wave function) of the condensed system; its use 
in studying vortex motion was pioneered by Fetter 0. 
The great advantage of the GP equation, as opposed to 
the hydrodynamic equations of a perfect fluid, is that it 
provides a detailed model for the vortex core, and may 
therefore be employed in the regime in which the vortex 
core size is comparable to or greater than other length 
scales in the problem. 



We elucidate, via this approach, how the velocity of 
the center of a vortex and the local fluid velocity dif- 
fer when the density of the condensate varies sufficiently 
rapidly in the neighborhood of a vortex - the situation in 
rapidly rotating condensates. Furthermore, we show how 
the difference between the vortex velocity and the fluid 
flow in such cases is compatible with Kelvin's theorem, a 
result with more general applications in hydrodynamics. 



THE VELOCITY OF A VORTEX 

The difference between the local fluid velocity and the 
velocity of the center of the vortex is brought out very 
clearly by the simple example of a single off-center vortex 
at radius & in a two-dimensional harmonic trap, V(r), 
of frequency u, in the limit in which the interparticle 
interaction is negligible. Such a system is described by 
a wave function, VP, that is a linear superposition of the 
(s-wave) oscillator ground state and a p-wave eigenstate 
of the trapping potential. At time t = 0, W ~ 
(C — &)e 1^1 l 2d , in the usual complex notation in which 
£ = x + iy; here d = ^Jh/muj is the oscillator length, and 
m the particle mass. Since the oscillator ground state 
has energy hto and the p-state has energy 2tku, the time 
dependent wave function is 

¥(£CV) ~ (e- 2luJt C~e- luJt b)e-^ 2/2d2 . (1) 

At time t the vortex is located at £ = be ZUJt , and thus 
its center precesses in the positive sense about the ori- 
gin at the trap frequency. On the other hand, the fluid 
velocity, given by (7i/m)V</>, where <f> is the phase of 
equals h/mp around the vortex, where p is the distance 
from the vortex. There is no background fluid flow; the 
only flow is that due to the vortex itself, and yet the vor- 
tex precesses with frequency oj. A naive application of 
Kelvin's theorem would suggest that the vortex should 
be stationary. 

We turn now to a more general analysis of the mo- 
tion of singly quantized vortices in two dimensions, as 
described by the time-dependent GP equation, 



ihd^/dt = -(h 2 /2m)V 2 ^ + + 52|*| 2 *, 



(2) 



where V(r ) is the transverse trapping potential, and g2 
is the two dimensional coupling constant. We normal- 
ize & by J d 2 r\^ 2 \ — N, where N is the total number 
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of particles. For a system that is uniform in the axial 
direction, 172 = 47r?i 2 a s /mZ, where a s is the s-wave in- 
teratomic scattering length, and Z is the height of the 
system. Given <!/, the left side, evaluated at the position 
of the vortex tells us the velocity of the vortex position. 
Thus the instantaneous velocity of the vortex center de- 
pends only on the value of V 2x P at the vortex position. 
The simplest case is that of a singly quantized cylindri- 
cally symmetric vortex at position Q in a spatially uni- 
form system, which is described by a wave function ^(C), 
~ C — d close to the vortex. More generally, the wave 
function of an asymmetric vortex, e.g., in an elliptic con- 
tainer, can include a term ~ (( — Q)* as well. Thus, to 
evaluate the GP equation at the vortex position, we can 
write the wave function in the neighborhood of a singly 
quantized vortex at position without loss of generality, 
as 



*(C, CV) = {[C -&(*)] + <* [C- C*Wn 



l+i4 



(3) 



where |a| < 1 for a vortex with positive circulation, and 
Q(C,C*^) an d < / , &(CiC*;0 are rea lj smooth functions at 
C = The background phase, <pb, is that of the su- 

perfluid with the singular contribution from the vortex 
at removed, while the function Q describes the back- 
ground density variation at the vortex. The background 
fluid velocity is Vb = Vb x + ivb y = 2(h,/m)d(f>b/d(*. Then 
the GP equation implies 



9 r -l. d r 
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■a-j£)(4> b + iQ)\ Ci . (4) 



This equation shows how the instantaneous velocity of 
the vortex center depends only on the slope and cur- 
vature of the wave function in the neighborhood of the 
vortex core; the velocity depends on the trapping poten- 
tial, V, and the interaction strength, g2, only insofar as 
they affect the value of the gradient of Q at the vortex 
position. 

After straightforward algebra we find that the veloc- 
ity, Ui, of the center of the vortex at Q is, in Cartesian 
coordinates, 



"Mix — Vbx 



1 t n -u9Q 



mRe(A x ) \ dy 



h fdQ .dQ 
u ly = v by H ^-r- I — 1mA 



mReA V d 



dy 



(5) 



where A = (1 — + a), and all quantities on the right 
are evaluated at Xi,yi. For a = 0, 



Ui = I Vb H z x VQ 

m 



(6) 



The vortex velocity is the sum of the local fluid velocity, 
Vb-, plus a correction, ~ VQ, of order (?i/m)Vrib, where 
rib ~ e _2< ^ describes the background density in the vortex 



core region. This correction is non-negligible if the back- 
ground density varies significantly over the vortex core. 
In the simple example above, v b = and the entire vortex 
velocity arises from the density gradient. As this calcula- 
tion illustrates, vortices do not in general move with the 
local fluid velocity, generated, e.g., by other vortices in 
the system. For the simple example in Eq. QJ, cf>b = 0, 
while VQ = r/d 2 , and Ui = ui(z x r.;). We note that 
these results hold for generalizations of the GP equation 
in which the interaction energy per particle is given by a 
local function of the density, as in hydrodynamics. 



HOW KELVIN'S THEOREM IS SATISFIED 

How can Eq. (0 be reconciled with Kelvin's theorem, 
that the circulation, § c v-ds around a contour, C, moving 
with the local fluid velocity, v(r), is conserved in time? 
(Here ds is the line element.) Kelvin's theorem applies 
to the GP equation, as well as the Schrodinger equation, 
since the flow described by these equations is potential. 
(The proof of Kelvin's theorem for the GP equation is 
given by Damski and Sacha @. As stressed in Ref. @ 
and earlier in Ref. the conditions for the proof of the 
theorem do not apply should the velocity be singular on 
the contour. This situation does not arise in the cases 
we consider here.) Since the theorem implies that a vor- 
tex cannot escape from a co-moving contour, no matter 
how small, the vortex would appear to be constrained 
to move with the background local fluid velocity. How- 
ever, as is illustrated in Fig. 1 for an off-center vortex 
in the harmonic oscillator potential, co-moving contours 
do not necessarily remain regular in time, but rather, 
due to the difference between the fluid motion about the 
vortex and the motion of the vortex center can become 
highly distorted. While in this example, the flow veloc- 
ity calculated from is always circular about the vortex 
center, the density gradient term in the vortex velocity 
which makes the vortex move with respect to the 
background flow, causes a mismatch between the motion 
of the contour and the motion \v\ — btu of the vortex 
center. An initially circular contour around the vortex 
begins to develop a lobe towards larger radii from the cen- 
ter. With time the lobe becomes extended and begins to 
wrap around. In Fig. 1 we show the fate, calculated nu- 
merically, of an initially circular contour centered on the 
vortex, at successive times within a quarter period, n/2u>. 
The standard argument that the vortex moves with the 
local fluid velocity assumes that the average velocity on 
the contour equals the fluid velocity at the center of the 
vortex. As this calculation indicates, it is not possible in 
general to deduce the motion of the vortex from Kelvin's 
theorem when the average velocity on the contour is no 
longer simply related to the vortex motion. 
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FIG. 1: Contour co-moving with the fluid, surrounding a sin- 
gle off-center vortex in a harmonic trap. The contour is shown 
in the stationary lab frame at successive equally spaced times, 
(a)-(d), within one-quarter of an oscillation period, n/2ui. 
Distances are measured in units of the oscillator length, d. 
The vortex is at radius d and the initial radius of the contour 
is 0.3d. Note the sharp fold of the contour in pane (d). 



APPLICATIONS 

We now apply the result Eq. (JSJ) to discuss more gen- 
eral vortex motion in detail in two limiting cases. The 
first is when the density scale height R is large com- 
pared with the vortex core size, determined in this case 
by the healing length, £; such a situation arises, e.g., in 
a system slowly rotating at angular frequency f2, when 
the interaction energy, ~ g^n (where n is the mean two- 
dimensional density), dominates the rotational energy, 
Ti£l, a strongly interacting regime in this sense. The 
second is the weakly-interacting rapidly-rotating lowest 
Landau level (LLL) limit in a harmonic trapping poten- 
tial, where g2n <C hQ. In the former case a vortex does 
move with the local fluid velocity, but as the harmonic 
oscillator example indicates, it does not in the latter case. 
The motion of vortices in the limit of large g2n may be ex- 
plicitly calculated by the method of asymptotic matching 
|l(lllllll2Lll3Lll4j . which we recall here. In this approach 
one solves the GP equation near the vortex by expanding 
the wave function about the solution for a uniform sys- 
tem to first order in the gradient of the potential. In the 
far field, at distances large compared with £, the veloc- 
ity is determined by solving the continuity equation with 
the density taken to be that in the absence of the vortex. 
Matching of the inner and outer solutions in the region 
where the radial coordinate p, measured with respect to 



the position of the vortex, is in the range ( « p « Ji, 
gives the vortex velocity and the parameters of the core 
wave function. The solution near the vortex core, which 
enters Eq. ©, has the form [3 El El El 



* = (l*o(p)|+x(p)cos0)e ie +"^ 



(7) 



where V^o is the vortex solution in a uniform medium, 
r\ and x are rea l functions of p, and the azimuthal an- 
gle 8, about the center of the vortex, is zero in the out- 
ward direction (opposite to the background density gra- 
dient near the vortex). Comparing this structure with 
the general form ©, we see that a vanishes in this limit 
and the corrections due to asymmetry of the vortex core 
are negligible; the background phase is <t>b = rj(p) sin 9. 
Using the analysis of Refs. El El El E3- we nn d 
X ~ p A at small p, and r\ ~ p. In addition, for small 
p, Q = — xcos^/l^ol ~ p 3 cos6. Thus, as we see from 
Eq. ©, the only contribution to the vortex velocity 
comes from the background phase of the wave function; 
the vortex moves with the background flow velocity. The 
back gro und fluid velocity, found from the far-field solu- 
tion fll E3I Eil Ei| , is given to logarithmic accuracy by 



MP) = ^-(2xVlnn)ln (■?-) , 
2m \R/ 



(8) 



in the region while in the core region, p/R should 

be replaced by £/R. Here n is the smooth (Thomas- 
Fermi) density outside the vortex core, given by Vn = 
— W(r)/g2- The flow velocity is induced by the vortex 
as a consequence of the density gradient in the system. 

In the rapidly rotating limit, on the other hand, in 
which the size of the vortex core is comparable to other 
lengths in the problem, vortices do not move with the 
background flow. In a harmonic trap the wave function 
in this limit is, to a first approximation, a linear superpo- 
sition of lowest Landau levels (LLL) in the Coriolis force 

El El, * ~ rij(C - 0)e HC|2/2d2 - Thus a = 0, and the 
velocity of the individual vortices, Eq. is given by, 



- Ci = *— MCI } 



(9) 



i.e., Ui = ujz x t\. In other words, each vortex precesses 
at the trap frequency, independent of the flow from other 
vortices! This result implies immediately that in any mo- 
tion of the vortices other than precession at the trap fre- 
quency, as for example in Tkachenko modes in the rapidly 
rotating limit 0] , higher Landau levels play a crucial role 
in the wave function. 

As a first step in describing the effects of interactions 
on vortex motion we consider the problem of steady state 
precession, in which there exists a frame rotating at an- 
gular velocity Q in which the system is stationary. The 
problem then is to determine f2 as a function of the in- 
teraction strength. In general, fl = d(H)/d(L), where 
L is the z component of the angular momentum. In 
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a system composed primarily of lowest Landau levels, 
(H) = uj((L) + N) + (H int ), so that Q-lu = d{H int )/d(L). 
As the angular momentum of the system increases, the 
particles spread out, lowering the average density; thus 
the derivative of the interaction energy with respect to 
the angular momentum is negative, implying that SI < uj. 
The first deviation of SI from uj is of order g 2 ■ 

Instead of calculating SI from d{H)/d(L), we discuss 
the precession of a single vortex (an effect measured, in 
the small core regime, by Anderson et al. 0]), directly in 
terms of Eq. (J5J) ■ Explicitly, we calculate the exact wave 
function of a single vortex to order b and g 2 , and from 
this show how the known vortex precession rate emerges 
for a vortex at small distance, 6, from the origin, to first 
order in the coupling strength |2fjj. (We set h=l in this 
section). In the frame rotating at SI, the GP equation 
assumes the form, 



(H -Ql-^ = -g 2 \-i>\H, 



(10) 



where Hq is the two-dimensional oscillator Hamiltonian, 
I — {Qd/dQ — (*d/d(*) is the angular momentum opera- 
tor, and /j, is the chemical potential in the rotating frame 
(equal to uj in the absence of interactions). Since SI — uj 
is first order in 52, as is fj, — uj, it follows from Eq. i|10|) 
that the interaction term mixes LLL components of or- 
der g 2 in \I/ for all angular momentum v. In the sense of 
quantum-mechanical perturbation theory, while the in- 
teraction is of order g 2 , the splitting of the LLL and thus 
energy denominators are also of order g 2 , leading to g® 
corrections. 

To calculate the vortex velocity from Eq. Q , or equiv- 
alently, (pi, one must also include the higher Landau level 
contributions to ^ explicitly; in general the exact wave 
function describing a vortex at position ( in terms of 
Landau levels is: 

*co = £ c va [p w (c, n - p*«,(Co, a)] e ~ |c|2/2d2 , 

ua 

(11) 

where the normalized Landau level wave functions are 
XW(C) = P w (C,C*)e HC|2/2d2 , with P va a polynomial, 
v = 0, ±1, ±2, ... is the angular momentum, and a — 
0,1,2,... is the radial quantum number. The state 
Xf<r(C) is an eigenstate of H with energy (2cr+ \v\ + 1)uj. 

The GP equation for a wave function stationary in the 
rotating frame, evaluated at the vortex position £ = b, is 



if 



to d(d( 

With given by Eq. (|llf> . this equation reduces to 



(12) 
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(n-uj^c^c—p^ 



dc 
(n + uj)c 



d 



d 



d(* to dCd(* 



P, 



(13) 



with all terms evaluated at b. The uj terms arise from 
de~'^' l 2d /dQdC,*. Since SI — uj is of order g 2 we need 
keep only LLL terms, er = 1 in the sum on the left, and 
of these only the v — 1 term is of order b°. On the right 
side, only higher Landau levels enter, since the P„o are 
independent of £*. 

The degree of mixing of higher v in the wave function 
of a single vortex depends on its distance b from the ori- 
gin. To leading order in b, C v q ~ for v > 1, while 
Coo ~ b. This relation is readily proven recursively: the 
interaction mixes terms \v0, Xw, and Xoo to give a v + 1 
component; since Coo begins at order b, then if C„o is 
of order W^ 1 , C + i.o begins at order b"; other mixings 
that give + 1 component lead to terms at least of or- 
der b v . An immediate consequence of this mixing is that 
even if one starts in the absence of interactions with a 
single vortex at b, the wave function in the presence of 
interactions describes further vortices at large distances, 
~ d 2 /b; in general, the LLL part of \& is proportional 
to a polynomial in £ of the form C,f(bC,/d 2 ) — b, to lead- 
ing order in b, where /(0) = 1, so that in addition to 
the zero at £ = 6, the polynomial has further zeroes at 
C| ~ d 2 jb. However, the particle density is negligible at 
such distances. 

As an explicit illustraion, let us start with a single vor- 
tex at position b described in the absence of interactions 
by 



*(C) 



N ( -b 



-|CI/2d 



(14) 



then with interactions the LLL wave function - even 
to order g 2 - includes levels (v, a — 0) with amplitude 
~ b v ~ x for v > 1. To first order in b, we must include 
as well the v — 2, a — LLL. Furthermore, the interac- 
tions to first order in g 2 mix in higher Landau levels with 
v = 0, 1, 2. However the v — 2 term contributes to the 
equation for the motion of the vortex only in order 6 2 , 
and thus can be neglected. Computing the coefficients 
Co- by expansion of the GP equation (|10|l one readily 
finds the exact wave function of a single vortex at posi- 
tion b, to first order in b and g 2 : 



mg 2 N 
IQtzTi 2 



(ICAf-2) 



+ •••} 



,-icr/2d 2 



(15) 



where C ~ y/ N/w is a normalization constant. The term 
~ C(IC/^| 2 — 2) is the v = l,a = 1 contribution. (Wc 
include the higher order b 2 term only to indicate that 
^(C) has a zero at £ = b.) 

Noting that a = here, we substitute this wave func- 
tion into Eq. (or equivalently, into Eq. ©) to derive 
the precession rate to first order in g 2 and lowest order 
in 6, El: 



SI — uj 



Ng 2 
'8nhd 2 ' 



(16) 
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This result illustrates the role of higher Landau levels 
in producing deviations of the precession frequency from 
the trap frequency. While the background velocity, v b , 
vanishes in the non-interacting case, we see in the pres- 
ence of interactions, from Eq. ||BJ|, that v b = bcoy/2, for a 
vortex centered at (x = b,y = 0), while the second term 
in JfjJ) equals (1 — mgN / A-Kh 2 )bujy /2. The precessional 
motion arises in part from the background flow and in 
part from the density gradient term. 

In addition, the chemical potential shift is given by 
fi — Uu) = 3<72^V/87T(i 2 , to first order in 52 and lowest 
order in b. To derive the b 2 corrections to the precession 
rate and chemical potential would require including the 
v = 3 components in <J . 



SUMMARY 

In summary, we have derived a general expression for 
the velocity of the center of a vortex in an atomic Bosc- 
Einstein condensate in terms of the background fluid ve- 
locity and the density gradients at the vortex center. The 
differences of the two velocities becomes important when 
the vortex core size becomes comparable to the scale 
height of the density variation. In particular, in the case 
of a rapidly rotating condensate in a harmonic trap, vor- 
tices precess at the transverse trap frequency, lo, minus 
corrections of order Na s /Z. 
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ments that led us to use the GP equation directly to 
study the motion of vortices, to Benny Lautrup for valu- 
able discussion, and Sandy Fetter and Bogdan Damski 
for useful correspondence. This research was supported 
in part by NSF Grants PHY03-55014 and PHY05-00914. 



[2 

[3 

[5. 
[6 

ir 

[8 

[9 

[10; 

[11 

[12 

[13 

[14 
[15 
[16 

[17. 
[18 

[19 

[20 



Schweikhard, V., Coddington, I., Engels, P., Mogen- 
dorff, V.P., & Cornell, E.A., (2004) Phys. Rev. Lett. 92, 
040404-1-4. 

Sinova, J., Hanna, C.B., & MacDonald, A.H. (2002) 
Phys. Rev. Lett. 89, 030403-1-4. 
Baym, G. (2003) Phys. Rev. Lett. 91, 110402-1-4. 
Thomson, W. (Lord Kelvin) (1867) Proc. Roy. Soc. (Ed- 
inburgh), 6, 94-105. 

Landau, L.D. & and Lifshitz, E.M. (1959) Fluid Mechan- 
ics (Addison- Wesley, Reading, Mass.). 
Fetter, AL. (1966) Phys. Rev. 151, 100-104. 
Damski, B. & Sacha, K. (2003) J. Phys. A: Math. Gen. 
36, 2339-2345. 

Koplik, J. & Levine, H. (1993) Phys. Rev. Lett. 71, 1375- 
1378. 

Rubinstein, BY. & Pismen, L.M. (1994) Physica D 78, 
1-10. 

Pismen, L.M. (1999) Vortices in Nonlinear Fields, (Ox- 
ford Univ. Press, Oxford). 

Svidzinsky, A. A. & Fetter, AL. (2000) Phys. Rev. A. 62, 
063617-1-14. 

Fetter, AL. & Svidzinsky, A. A. (2001) J. Phys.: Con- 
dens. Matt. 13, R135-R194. 

Anglin, J.R. (2002) Phys. Rev. A 65, 063611-1-14. 
Lundh, E. & Ao, P. (2000) Phys. Rev. A 61, 063612-1-7. 
Sheehy, D.E. & Radzihovsky, L. (2004) Phys. Rev. A 70, 
051602(R)-l-4. 

Ho, T.-L. (2001) Phys. Rev. Lett. 87, 060403-1-4. 
Watanabe, G., Baym, G, & Pethick, C.J. (2004) Phys. 
Rev. Lett. 93, 190401-1-4. 

Anderson, B.P., Haljan, P.C., Wieman, C.E., & Cornell, 

E.A. (2000) Phys. Rev. Lett. 85, 2857-2860. 

Linn, M. & Fetter, AL. (2000) Phys. Rev. A 61, 063603- 

1-4. 



[1] Coddington, I., Engels, P., Schweikhard, V., & Cornell, 
E.A. (2003) Phys. Rev. Lett. 91, 100402-1-4. 



